gtex_CT_TS_TW_h2
Heather E. Wheeler
2015-10-20 16:13:30
GTEx OTD joint heritability. Local h2 is estimated with SNPs within 1 Mb of each gene. distal h2 is estimated with SNPs that are eQTLs in the Framingham Heart Study on other chromosomes (FDR < 0.05).
library(ggplot2)
library(reshape2)
library(dplyr)
library(tidyr)
library(GGally)
library(grid)
"%&%" = function(a,b) paste(a,b,sep="")
source('/Volumes/im-lab/nas40t2/hwheeler/PrediXcan-Paper/scripts/Heather/make-figures/multiplot.R')
my.dir <- '/Volumes/im-lab/nas40t2/hwheeler/cross-tissue/'
fig.dir <- '~/GitHub/GenArch/GenArchPaper/Figures/'Tissue-Specific joint local estimates
tislist <- scan(my.dir %&% 'TS.ten.tissue.list',sep="\n",what="character")
ts <- data.frame()
set.seed(12345)
for(tis in tislist){
data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")
#replace NA with h2=0 and se=mean(se) #change h2 to mean h2?
loc.jt <- select(data,tissue,loc.jt.h2,loc.jt.se) %>% mutate(loc.jt.h2=ifelse(is.na(loc.jt.h2),0,loc.jt.h2), loc.jt.se=ifelse(is.na(loc.jt.se),sample(loc.jt.se[is.na(loc.jt.se)==FALSE],size=length(loc.jt.se[is.na(loc.jt.se)==TRUE])),loc.jt.se)) %>% mutate(ymin = pmax(0, loc.jt.h2 - 2 * loc.jt.se), ymax = pmin(1, loc.jt.h2 + 2 * loc.jt.se)) %>% arrange(loc.jt.h2) %>% mutate(nonzeroCI=ymin>0,place=1:length(data$tissue))
ts <- rbind(ts,loc.jt)
}
p<-ggplot(ts,aes(x=place,y=loc.jt.h2,ymin=ymin,ymax=ymax,color=nonzeroCI) ) + facet_wrap(~tissue,ncol=2) + geom_pointrange(col='gray')+geom_point()+ ylab(expression("local h"^2)) + xlab(expression("genes ordered by joint local h"^2))+theme_bw() + coord_cartesian(ylim=c(0,1)) + ggtitle("CT and TS")
###calc % nonzero for each tissue
### ADD to plot legend
pvec<-vector()
h2vec<-vector()
a<-ts %>% select(tissue,nonzeroCI) %>% spread(tissue,nonzeroCI)
for(i in 1:10){
tis<-colnames(a)[i]
cat("\n\n---",tis,"---\n")
print(table(a[,i]))
per <- signif(table(a[,i])/sum(table(a[,i])),3)*100
print(per)
pvec <- c(pvec,per[2])
}##
##
## --- cross-tissue ---
##
## FALSE TRUE
## 14937 2014
##
## FALSE TRUE
## 88.1 11.9
##
##
## --- Adipose - Subcutaneous ---
##
## FALSE TRUE
## 16851 100
##
## FALSE TRUE
## 99.40 0.59
##
##
## --- Artery - Tibial ---
##
## FALSE TRUE
## 16808 143
##
## FALSE TRUE
## 99.200 0.844
##
##
## --- Heart - Left Ventricle ---
##
## FALSE TRUE
## 16827 124
##
## FALSE TRUE
## 99.300 0.732
##
##
## --- Lung ---
##
## FALSE TRUE
## 16853 98
##
## FALSE TRUE
## 99.400 0.578
##
##
## --- Muscle - Skeletal ---
##
## FALSE TRUE
## 16740 211
##
## FALSE TRUE
## 98.80 1.24
##
##
## --- Nerve - Tibial ---
##
## FALSE TRUE
## 16767 184
##
## FALSE TRUE
## 98.90 1.09
##
##
## --- Skin - Sun Exposed (Lower leg) ---
##
## FALSE
## 16951
##
## FALSE
## 100
##
##
## --- Thyroid ---
##
## FALSE
## 16951
##
## FALSE
## 100
##
##
## --- Whole Blood ---
##
## FALSE TRUE
## 16707 244
##
## FALSE TRUE
## 98.60 1.44###calc mean h2 for each tissue
a<-ts %>% select(tissue,loc.jt.h2) %>% spread(tissue,loc.jt.h2)
for(i in 1:10){
tis<-colnames(a)[i]
meanh2 <- signif(mean(a[,i],na.rm=TRUE),3)
cat(tis,"mean h2:",meanh2,"\n")
h2vec <- c(h2vec,meanh2)
}## cross-tissue mean h2: 0.0406
## Adipose - Subcutaneous mean h2: 0.0128
## Artery - Tibial mean h2: 0.015
## Heart - Left Ventricle mean h2: 0.0181
## Lung mean h2: 0.0134
## Muscle - Skeletal mean h2: 0.0135
## Nerve - Tibial mean h2: 0.017
## Skin - Sun Exposed (Lower leg) mean h2: 0
## Thyroid mean h2: 0.0514
## Whole Blood mean h2: 0.1pvec<-ifelse(is.na(pvec),0,pvec)
ann_text <- data.frame( loc.jt.h2 = rep(0.78,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2<-p+geom_text(data=ann_text,aes(label=paste("mean_h^2 ==",mean_h2,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3)
ann_text <- data.frame( loc.jt.h2 = rep(0.9,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2+geom_text(data=ann_text,aes(label=paste("percent_TRUE ==",percent,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3)+ theme(legend.justification=c(0,1), legend.position=c(0,1))
Tissue-Specific joint global estimates
tislist <- scan(my.dir %&% 'TS.ten.tissue.list',sep="\n",what="character")
ts <- data.frame()
for(tis in tislist){
data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")
glo.jt <- select(data,tissue,glo.jt.h2,glo.jt.se) %>% mutate(glo.jt.h2=ifelse(is.na(glo.jt.h2),0,glo.jt.h2), glo.jt.se=ifelse(is.na(glo.jt.se),sample(glo.jt.se[is.na(glo.jt.se)==FALSE],size=length(glo.jt.se[is.na(glo.jt.se)==TRUE])),glo.jt.se)) %>% mutate(ymin = pmax(0, glo.jt.h2 - 2 * glo.jt.se), ymax = pmin(1, glo.jt.h2 + 2 * glo.jt.se)) %>% arrange(glo.jt.h2) %>% mutate(nonzeroCI=ymin>0,place=1:length(data$tissue))
ts <- rbind(ts,glo.jt)
}
p<-ggplot(ts,aes(x=place,y=glo.jt.h2,ymin=ymin,ymax=ymax,color=nonzeroCI) ) + facet_wrap(~tissue,ncol=2) + geom_pointrange(col='gray')+geom_point()+ylab(expression("global h"^2)) + xlab(expression("genes ordered by joint global h"^2))+theme_bw() + coord_cartesian(ylim=c(0,1)) + ggtitle("CT and TS")
###calc % nonzero for each tissue
### ADD to plot legend
pvec<-vector()
h2vec<-vector()
a<-ts %>% select(tissue,nonzeroCI) %>% spread(tissue,nonzeroCI)
for(i in 1:10){
tis<-colnames(a)[i]
cat("\n\n---",tis,"---\n")
print(table(a[,i]))
per <- signif(table(a[,i])/sum(table(a[,i])),3)*100
print(per)
pvec <- c(pvec,per[2])
}##
##
## --- cross-tissue ---
##
## FALSE TRUE
## 16618 333
##
## FALSE TRUE
## 98.00 1.96
##
##
## --- Adipose - Subcutaneous ---
##
## FALSE TRUE
## 16914 37
##
## FALSE TRUE
## 99.800 0.218
##
##
## --- Artery - Tibial ---
##
## FALSE TRUE
## 16937 14
##
## FALSE TRUE
## 99.9000 0.0826
##
##
## --- Heart - Left Ventricle ---
##
## FALSE
## 16951
##
## FALSE
## 100
##
##
## --- Lung ---
##
## FALSE TRUE
## 16944 7
##
## FALSE TRUE
## 100.0000 0.0413
##
##
## --- Muscle - Skeletal ---
##
## FALSE TRUE
## 16764 187
##
## FALSE TRUE
## 98.9 1.1
##
##
## --- Nerve - Tibial ---
##
## FALSE TRUE
## 16950 1
##
## FALSE TRUE
## 1.0e+02 5.9e-03
##
##
## --- Skin - Sun Exposed (Lower leg) ---
##
## FALSE TRUE
## 16836 115
##
## FALSE TRUE
## 99.300 0.678
##
##
## --- Thyroid ---
##
## FALSE TRUE
## 16939 12
##
## FALSE TRUE
## 99.9000 0.0708
##
##
## --- Whole Blood ---
##
## FALSE TRUE
## 16825 126
##
## FALSE TRUE
## 99.300 0.743###calc mean h2 for each tissue
a<-ts %>% select(tissue,glo.jt.h2) %>% spread(tissue,glo.jt.h2)
for(i in 1:10){
tis<-colnames(a)[i]
meanh2 <- signif(mean(a[,i],na.rm=TRUE),3)
cat(tis,"mean h2:",meanh2,"\n")
h2vec <- c(h2vec,meanh2)
}## cross-tissue mean h2: 0.107
## Adipose - Subcutaneous mean h2: 0.092
## Artery - Tibial mean h2: 0.117
## Heart - Left Ventricle mean h2: 0.162
## Lung mean h2: 0.107
## Muscle - Skeletal mean h2: 0.0895
## Nerve - Tibial mean h2: 0.149
## Skin - Sun Exposed (Lower leg) mean h2: 0.151
## Thyroid mean h2: 0.11
## Whole Blood mean h2: 0.108pvec<-ifelse(is.na(pvec),0,pvec)
ann_text <- data.frame( glo.jt.h2 = rep(0.78,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2<-p+geom_text(data=ann_text,aes(label=paste("mean_h^2 ==",mean_h2,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3)
ann_text <- data.frame( glo.jt.h2 = rep(0.9,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2+geom_text(data=ann_text,aes(label=paste("percent_TRUE ==",percent,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3)+ theme(legend.justification=c(0,1), legend.position=c(0,1))
Tissue-Wide joint local estimates
tislist <- scan(my.dir %&% 'TW.ten.tissue.list',sep="\n",what="character")
ts <- data.frame()
for(tis in tislist){
data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")
loc.jt <- select(data,tissue,loc.jt.h2,loc.jt.se) %>% mutate(loc.jt.h2=ifelse(is.na(loc.jt.h2),0,loc.jt.h2), loc.jt.se=ifelse(is.na(loc.jt.se),sample(loc.jt.se[is.na(loc.jt.se)==FALSE],replace=TRUE,size=length(loc.jt.se[is.na(loc.jt.se)==TRUE])),loc.jt.se)) %>% mutate(ymin = pmax(0, loc.jt.h2 - 2 * loc.jt.se), ymax = pmin(1, loc.jt.h2 + 2 * loc.jt.se)) %>% arrange(loc.jt.h2) %>% mutate(nonzeroCI=ymin>0,place=1:length(data$tissue))
ts <- rbind(ts,loc.jt)
}
p<-ggplot(ts,aes(x=place,y=loc.jt.h2,ymin=ymin,ymax=ymax,color=nonzeroCI) ) + facet_wrap(~tissue,ncol=2) + geom_pointrange(col='gray')+geom_point()+ylab(expression("local h"^2)) + xlab(expression("genes ordered by joint local h"^2))+theme_bw() + coord_cartesian(ylim=c(0,1)) + ggtitle("CT and TW")
###calc % nonzero for each tissue
### ADD to plot legend
pvec<-vector()
h2vec<-vector()
a<-ts %>% select(tissue,nonzeroCI) %>% spread(tissue,nonzeroCI)
for(i in 1:10){
tis<-colnames(a)[i]
cat("\n\n---",tis,"---\n")
print(table(a[,i]))
per <- signif(table(a[,i])/sum(table(a[,i])),3)*100
print(per)
pvec <- c(pvec,per[2])
}##
##
## --- cross-tissue ---
##
## FALSE TRUE
## 14937 2014
##
## FALSE TRUE
## 88.1 11.9
##
##
## --- Adipose-Subcutaneous ---
##
## FALSE TRUE
## 17882 588
##
## FALSE TRUE
## 96.80 3.18
##
##
## --- Artery-Tibial ---
##
## FALSE TRUE
## 17745 674
##
## FALSE TRUE
## 96.30 3.66
##
##
## --- Heart-LeftVentricle ---
##
## FALSE TRUE
## 17882 411
##
## FALSE TRUE
## 97.80 2.25
##
##
## --- Lung ---
##
## FALSE TRUE
## 17949 523
##
## FALSE TRUE
## 97.20 2.83
##
##
## --- Muscle-Skeletal ---
##
## FALSE TRUE
## 17893 606
##
## FALSE TRUE
## 96.70 3.28
##
##
## --- Nerve-Tibial ---
##
## FALSE TRUE
## 17691 800
##
## FALSE TRUE
## 95.70 4.33
##
##
## --- Skin-SunExposed(Lowerleg) ---
##
## FALSE TRUE
## 17862 622
##
## FALSE TRUE
## 96.60 3.37
##
##
## --- Thyroid ---
##
## FALSE TRUE
## 17746 724
##
## FALSE TRUE
## 96.10 3.92
##
##
## --- WholeBlood ---
##
## FALSE TRUE
## 17889 589
##
## FALSE TRUE
## 96.80 3.19###calc mean h2 for each tissue
a<-ts %>% select(tissue,loc.jt.h2) %>% spread(tissue,loc.jt.h2)
for(i in 1:10){
tis<-colnames(a)[i]
meanh2 <- signif(mean(a[,i],na.rm=TRUE),3)
cat(tis,"mean h2:",meanh2,"\n")
h2vec <- c(h2vec,meanh2)
}## cross-tissue mean h2: 0.0406
## Adipose-Subcutaneous mean h2: 0.0172
## Artery-Tibial mean h2: 0.0186
## Heart-LeftVentricle mean h2: 0.013
## Lung mean h2: 0.00396
## Muscle-Skeletal mean h2: 0.0157
## Nerve-Tibial mean h2: 0.0204
## Skin-SunExposed(Lowerleg) mean h2: 0.0172
## Thyroid mean h2: 0.0199
## WholeBlood mean h2: 0.015pvec<-ifelse(is.na(pvec),0,pvec)
ann_text <- data.frame( loc.jt.h2 = rep(0.78,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2<-p+geom_text(data=ann_text,aes(label=paste("mean_h^2 ==",mean_h2,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3)
ann_text <- data.frame( loc.jt.h2 = rep(0.9,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2+geom_text(data=ann_text,aes(label=paste("percent_TRUE ==",percent,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3)+ theme(legend.justification=c(0,1), legend.position=c(0,1))
Tissue-Wide joint global estimates
tislist <- scan(my.dir %&% 'TW.ten.tissue.list',sep="\n",what="character")
ts <- data.frame()
for(tis in tislist){
data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")
glo.jt <- select(data,tissue,glo.jt.h2,glo.jt.se) %>% mutate(glo.jt.h2=ifelse(is.na(glo.jt.h2),0,glo.jt.h2), glo.jt.se=ifelse(is.na(glo.jt.se),sample(glo.jt.se[is.na(glo.jt.se)==FALSE],replace=TRUE,size=length(glo.jt.se[is.na(glo.jt.se)==TRUE])),glo.jt.se)) %>% mutate(ymin = pmax(0, glo.jt.h2 - 2 * glo.jt.se), ymax = pmin(1, glo.jt.h2 + 2 * glo.jt.se)) %>% arrange(glo.jt.h2) %>% mutate(nonzeroCI=ymin>0,place=1:length(data$tissue))
ts <- rbind(ts,glo.jt)
}
p<-ggplot(ts,aes(x=place,y=glo.jt.h2,ymin=ymin,ymax=ymax,color=nonzeroCI) ) + facet_wrap(~tissue,ncol=2) + geom_pointrange(col='gray')+geom_point()+ylab(expression("global h"^2)) + xlab(expression("genes ordered by joint global h"^2))+theme_bw() + coord_cartesian(ylim=c(0,1)) + ggtitle("CT and TW")
###calc % nonzero for each tissue
### ADD to plot legend
pvec<-vector()
h2vec<-vector()
a<-ts %>% select(tissue,nonzeroCI) %>% spread(tissue,nonzeroCI)
for(i in 1:10){
tis<-colnames(a)[i]
cat("\n\n---",tis,"---\n")
print(table(a[,i]))
per <- signif(table(a[,i])/sum(table(a[,i])),3)*100
print(per)
pvec <- c(pvec,per[2])
}##
##
## --- cross-tissue ---
##
## FALSE TRUE
## 16618 333
##
## FALSE TRUE
## 98.00 1.96
##
##
## --- Adipose-Subcutaneous ---
##
## FALSE TRUE
## 18454 16
##
## FALSE TRUE
## 99.9000 0.0866
##
##
## --- Artery-Tibial ---
##
## FALSE TRUE
## 18415 4
##
## FALSE TRUE
## 100.0000 0.0217
##
##
## --- Heart-LeftVentricle ---
##
## FALSE
## 18293
##
## FALSE
## 100
##
##
## --- Lung ---
##
## FALSE TRUE
## 18469 3
##
## FALSE TRUE
## 100.0000 0.0162
##
##
## --- Muscle-Skeletal ---
##
## FALSE TRUE
## 18355 144
##
## FALSE TRUE
## 99.200 0.778
##
##
## --- Nerve-Tibial ---
##
## FALSE TRUE
## 18490 1
##
## FALSE TRUE
## 1.00e+02 5.41e-03
##
##
## --- Skin-SunExposed(Lowerleg) ---
##
## FALSE TRUE
## 18458 26
##
## FALSE TRUE
## 99.900 0.141
##
##
## --- Thyroid ---
##
## FALSE
## 18470
##
## FALSE
## 100
##
##
## --- WholeBlood ---
##
## FALSE TRUE
## 18325 153
##
## FALSE TRUE
## 99.200 0.828###calc mean h2 for each tissue
a<-ts %>% select(tissue,glo.jt.h2) %>% spread(tissue,glo.jt.h2)
for(i in 1:10){
tis<-colnames(a)[i]
meanh2 <- signif(mean(a[,i],na.rm=TRUE),3)
cat(tis,"mean h2:",meanh2,"\n")
h2vec <- c(h2vec,meanh2)
}## cross-tissue mean h2: 0.107
## Adipose-Subcutaneous mean h2: 0.0478
## Artery-Tibial mean h2: 0.051
## Heart-LeftVentricle mean h2: 0.0407
## Lung mean h2: 0.0729
## Muscle-Skeletal mean h2: 0.0421
## Nerve-Tibial mean h2: 0.0448
## Skin-SunExposed(Lowerleg) mean h2: 0.0428
## Thyroid mean h2: 0
## WholeBlood mean h2: 0.05pvec<-ifelse(is.na(pvec),0,pvec)
ann_text <- data.frame( glo.jt.h2 = rep(0.78,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2<-p+geom_text(data=ann_text,aes(label=paste("mean_h^2 ==",mean_h2,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3)
ann_text <- data.frame( glo.jt.h2 = rep(0.9,10), place = rep(5000,10), percent= pvec, mean_h2 = h2vec, tissue = factor(colnames(a)), ymin=rep(0.9,10), ymax=rep(0.9,10), nonzeroCI=rep(NA,10), se=rep(0.9,10))
p2+geom_text(data=ann_text,aes(label=paste("percent_TRUE ==",percent,sep="")),color="black",show_guide=F,parse=T,hjust=0,size=3)+ theme(legend.justification=c(0,1), legend.position=c(0,1))
Cross-tissue v. Tissue-specific joint h2 and se
tislist <- scan(my.dir %&% 'ten.tissue.list',sep="\n",what="character")
ts <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.' %&% tislist[1] %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")
##LOCAL
tsh2 <-ts %>% select(ensid,loc.jt.h2)
colnames(tsh2) = c("ensid",tislist[1])
tsse <-ts %>% select(ensid,loc.jt.se)
colnames(tsse) = c("ensid",tislist[1])
for(tis in tislist[2:10]){
data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.TS.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")
datah2 <-data %>% select(ensid,loc.jt.h2)
colnames(datah2) = c("ensid",tis)
datase <-data %>% select(ensid,loc.jt.se)
colnames(datase) = c("ensid",tis)
tsh2 <- inner_join(tsh2,datah2,by="ensid")
tsse <- inner_join(tsse,datase,by="ensid")
}
gtsh2<-gather(tsh2,"CrossTissue","Tissue",3:11)
colnames(gtsh2) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtsh2,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red') + ylab(expression('Cross-Tissue h'^2)) + xlab(expression('Tissue-Specific Joint local h'^2)) + coord_cartesian(ylim=c(-0.1,1.1),xlim=c(-0.1,1.1)) + theme_bw() 
gtsse<-gather(tsse,"CrossTissue","Tissue",3:11)
colnames(gtsse) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtsse,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red') + ylab('Cross-Tissue SE') + xlab('Tissue-Specific Joint local SE') + coord_cartesian(ylim=c(-0.01,0.18),xlim=c(-0.01,0.18)) + theme_bw() 
##GLOBAL
tsh2 <-ts %>% select(ensid,glo.jt.h2)
colnames(tsh2) = c("ensid",tislist[1])
tsse <-ts %>% select(ensid,glo.jt.se)
colnames(tsse) = c("ensid",tislist[1])
for(tis in tislist[2:10]){
data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.TS.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")
datah2 <-data %>% select(ensid,glo.jt.h2)
colnames(datah2) = c("ensid",tis)
datase <-data %>% select(ensid,glo.jt.se)
colnames(datase) = c("ensid",tis)
tsh2 <- inner_join(tsh2,datah2,by="ensid")
tsse <- inner_join(tsse,datase,by="ensid")
}
gtsh2<-gather(tsh2,"CrossTissue","Tissue",3:11)
colnames(gtsh2) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtsh2,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red') + ylab(expression('Cross-Tissue h'^2)) + xlab(expression('Tissue-Specific Joint global h'^2)) + coord_cartesian(ylim=c(-0.1,1.1),xlim=c(-0.1,1.1)) + theme_bw() 
gtsse<-gather(tsse,"CrossTissue","Tissue",3:11)
colnames(gtsse) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtsse,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red') + ylab('Cross-Tissue SE') + xlab('Tissue-Specific Joint global SE') + coord_cartesian(ylim=c(-0.01,1.4),xlim=c(-0.01,1.4)) + theme_bw() 
Cross-tissue v. Tissue-wide joint h2 and se
tislist <- scan(my.dir %&% 'rmTW.ten.tissue.list',sep="\n",what="character")
tw <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.' %&% tislist[1] %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")
##LOCAL
twh2 <-tw %>% select(ensid,loc.jt.h2)
colnames(twh2) = c("ensid",tislist[1])
twse <-tw %>% select(ensid,loc.jt.se)
colnames(twse) = c("ensid",tislist[1])
for(tis in tislist[2:10]){
data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.TW.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")
datah2 <-data %>% select(ensid,loc.jt.h2)
colnames(datah2) = c("ensid",tis)
datase <-data %>% select(ensid,loc.jt.se)
colnames(datase) = c("ensid",tis)
twh2 <- inner_join(twh2,datah2,by="ensid")
twse <- inner_join(twse,datase,by="ensid")
}
gtwh2<-gather(twh2,"CrossTissue","Tissue",3:11)
colnames(gtwh2) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtwh2,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red') + ylab(expression('Cross-Tissue h'^2)) + xlab(expression('Tissue-Wide Joint local h'^2)) + coord_cartesian(ylim=c(-0.1,1.1),xlim=c(-0.1,1.1)) + theme_bw() 
gtwse<-gather(twse,"CrossTissue","Tissue",3:11)
colnames(gtwse) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtwse,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red') + ylab('Cross-Tissue SE') + xlab('Tissue-Wide Joint local SE') + coord_cartesian(ylim=c(-0.01,0.21),xlim=c(-0.01,0.21)) + theme_bw() 
##GLOBAL
twh2 <-tw %>% select(ensid,glo.jt.h2)
colnames(twh2) = c("ensid",tislist[1])
twse <-tw %>% select(ensid,glo.jt.se)
colnames(twse) = c("ensid",tislist[1])
for(tis in tislist[2:10]){
data <- read.table(my.dir %&% 'gtex-h2-estimates/GTEx.TW.' %&% tis %&% '.h2.all.models_FHSfdr0.05.Chr1-22_globalOtherChr.2015-10-06.txt',header=T,sep="\t")
datah2 <-data %>% select(ensid,glo.jt.h2)
colnames(datah2) = c("ensid",tis)
datase <-data %>% select(ensid,glo.jt.se)
colnames(datase) = c("ensid",tis)
twh2 <- inner_join(twh2,datah2,by="ensid")
twse <- inner_join(twse,datase,by="ensid")
}
gtwh2<-gather(twh2,"CrossTissue","Tissue",3:11)
colnames(gtwh2) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtwh2,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red') + ylab(expression('Cross-Tissue h'^2)) + xlab(expression('Tissue-Wide Joint global h'^2)) + coord_cartesian(ylim=c(-0.1,1.1),xlim=c(-0.1,1.1)) + theme_bw() 
gtwse<-gather(twse,"CrossTissue","Tissue",3:11)
colnames(gtwse) <- c('ensid','CrossTissue','TissueName','Tissue')
ggplot(gtwse,aes(x=Tissue,y=CrossTissue)) +facet_wrap(~TissueName,scales="fixed",ncol=3) + geom_point(alpha=1/5) + geom_abline(intercept=0, slope=1,color='red') + ylab('Cross-Tissue SE') + xlab('Tissue-Wide Joint global SE') + coord_cartesian(ylim=c(-0.01,1.6),xlim=c(-0.01,1.6)) + theme_bw() 